Origin of Subgap States in Normal-Insulator-Superconductor van der Waals Heterostructures

Superconductivity in van der Waals materials, such as NbSe2 and TaS2, is fundamentally novel due to the effects of dimensionality, crystal symmetries, and strong spin–orbit coupling. In this work, we perform tunnel spectroscopy on NbSe2 by utilizing MoS2 or hexagonal boron nitride (hBN) as a tunnel barrier. We observe subgap excitations and probe their origin by studying various heterostructure designs. We show that the edge of NbSe2 hosts many defect states, which strongly couple to the superconductor and form Andreev bound states. Furthermore, by isolating the NbSe2 edge we show that the subgap states are ubiquitous in MoS2 tunnel barriers but absent in hBN tunnel barriers, suggesting defects in MoS2 as their origin. Their magnetic nature reveals a singlet- or a doublet-type ground state, and based on nearly vanishing g factors or avoided crossings of subgap excitations, we highlight the role of strong spin–orbit coupling.

all devices also contain a graphite crystal transferred over NbSe 2 , which acts as an ohmic contact to NbSe 2 . Normal contact regions are defined by ebeam lithography and Ti/Au (5 nm/50 nm) is deposited by ebeam evaporation to create the tunnel junctions. The area of the tunnel junctions is typically ∼ 1 µm 2 to 3 µm 2 , unless stated otherwise. All details of device parameters can be found in the device table at the end of the Supporting Information (SI).
The measurements were performed in a 3 He fridge at a base temperature of ∼ 250 mK. The electrical lines are filtered by using pi-filters at the breakout box and tapeworm filters at the cold finger. A small ac voltage (V ac < k B T , where k B is the Boltzmann constant and T is the cryostat temperature) is added to the biasing dc voltage V using a transformer and the lockin amplifier 2 records the output current via an external current-voltage amplifier. Due to a cryostat installation error the actual magnetic fields may be smaller by up to 10 % than reported here.

II. TUNNEL SPECTRA FITTING
In the main text we presented a fit of the tunnel spectra with simple BCS-type density of states, shown in Figure 1b also reproduced here in Figure S1a. However, the tunnel spectra of NbSe 2 is known to exhibit features that deviate from the simple BCS model (an isotropic, single gap). Other models have been employed to fit the tunnel spectra, and two important ones are 1) The two-band model [1,2], for which a fit is shown in Figure S1b and 2) the anisotropic-gap (single gap) model [3] with a fit in Figure S1c. The temperature for all fits is fixed to the measurement temperature, T = 255 mK. In our view the two band model presents the most consistent picture, where the superconductivity in thicker NbSe 2 is thought to arise from multiple bands. With the layer number decreasing in NbSe 2 the superconducting gap approaches a single gap model and is consistent with bandstructure calculations. In contrast it is hard to justify the change of anisotropy as the NbSe 2 becomes thinner.   Two band model with parameters ∆ 1 = 1.1 meV, Γ 1 = 0.08 meV, ∆ 2 = 0.9 meV, Γ 2 = 0.08 meV and N 1 /N 2 = 0.6.

III. TEMPERATURE DEPENDENCE OF THE TUNNEL SPECTRA
Here we highlight the role of superconductivity by showing temperature dependent tunnel spectra. The tunnel spectra exhibits broadening with temperature, as shown in the differential conductance map in Figure S2a. The superconducting gap is no longer visible for T ≥ 4.4 K and is reasonable for a thin (∼ 3 nm) NbSe 2 . In addition, the subgap excitations visible at a base temperature T = 270 mK are no longer visible at T = 1.5 K, as shown in Figure S2b. Similarly, for another device D21, shown in Figure S2c, the subgap excitations are no longer visible at a larger temperature T = 4 K, and the tunnel spectrum exhibits thermal broadening.  While measuring in the wide range of magnetic field ( Figure S4 and Figure S5), we take care of conductance jumps due to vortices by shifting dI/dV. The shift was done such that at all magnetic field values, the conductance at V b = 0 has the same value. Such a wide peak hinder us to distinguish if the lines cross as in Figure S4a or if the dI/dV peak sticks to V b =0 above 6 T.

VII. DEVICE DETAILS
The